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An alternative methodology to investigate indirect polyatomic processes with quasi-classical trajec- 
tories is proposed, which effectively avoids any binning or weighting procedure whiie provides rovi- 
brational resofution. Initiai classical states are started in terms of angle-action variables to closely 
match the quantum experimental conditions and later transformed into Cartesian coordinates, fol- 
lowing an algorithm very recently published [J. Chem. Phys. 130, 114103 (2009)]. Trajectories are 
then propagated using the 'association' picture, i.e. an inverse dynamics simulation in the spirit of 
the exit-channel corrected phase space theory of Hamilton and Brumer [J. Chem. Phys. 82, 595 
(1985)], which is shown to be particularly convenient. Finally, an approximate quasi-classical for- 
mula is provided which under general conditions can be used to add possible rotational structures 
into the vibrationally-resolved quasi-classical distributions. To introduce the method and illustrate 
its capabilities, correlated translational energy distributions from recent experiments in the photo- 
dissociation of ketene at 308 nm [J. Chem. Phys. 124, 014303 (2006)] are investigated. Quite 
generally, the overall theoretical algorithm reduces the total number of trajectories to integrate and 
allows for fully theoretical predictions of experiments on polyatomics. 

PACS numbers: 31.15.xg, 34.10.-l-x, 82.20.Bc, 82.37.Np 



I. INTRODUCTION 

The 20th century witnessed both the advent of modern 
computers and the associated exponential growth in com- 
putational power, which continues nowadays. As more 
can be done in significantly less time, simulation of in- 
herently more complex processes and/or larger systems is 
attempted. In the field of chemical physics, what started 
as describing atomic and diatomic molecular processes 
currently looks after the simulation of, e.g. larger poly- 
atomics, clusters, phenomena at surfaces and biomolecu- 
lar complexes of interest. 

However, despite all the effort being devoted to develop 
approaches which include quantum effects in the dynam- 
ics of large molecular systems realistic quantum- 
mechanical simulations of polyatomic processes remain 
prohibitively expensive. To date, the only general theo- 
retical alternative is the quasi-classical trajectory method 
(QCT) 0, H, 0|i which is based on numerically solving 
Hamilton's (or Newton's) equations of motion. Primar- 
ily, this is due to the cumbersome basis sets that need to 
be handled within quantum calculations in order to avoid 
crude approximations, but also partially because devel- 
oping general, i.e. scalable classical trajectory codes is 
relatively easier. 

To allow for comparison with experiments, both the 
selection of classical initial conditions and the statisti- 
cal treatment of QCT results should be semiclassical in 



* Corresponding author. lE-mail: mleo@instec.cu1 

t Corresponding author. lE-mail: l.bonnet@ism.u-bordeauxl.frl 



spirit. In particular, to mimic quantum features using 
classical (continuous) magnitudes, additional binning or 
weighting procedures are quite difficult to avoid. As 
implied, the underlying source of the problem — as with 
some other drawbacks of QCT — is that classically, vi- 
brations and rotations are not quantised per se. In 
this regard, important advances have been made in the 
last few years by replacing the standard/histogram bin- 
"^^"^K, ^^p'X wpcedvLVQ by a Gaussian weighting (GW) 
[1 i, E i, S [M [O, [l2- In the latter, the different 
trajectories are assigned different (Gaussian) statistical 
weights, larger the closer its final actions to integer val- 
ues. Despite its many promising capabilities, we have re- 
cently pointed out how the convergence of GW could be 
affected when applied to processes involving polyatomics 
prj . This is mainly because increasing the number of 
atoms in a molecule will almost certainly increase the 
number of modes for which a quantum-mechanical de- 
scription is convenient — and Gaussian weights must be 
assigned. As a result, the overall weight of a given trajec- 
tory could be significantly reduced, in particular if several 
actions vary simultaneously. 

At the same time, powerful experimental techniques 
have been developed during the precedent decades and 
very precise measurements of state-to-state observables 
are now possible. Correlated product distributions — 
where final state and/or translational energy distribu- 
tions are measured in coincidence with given rovibra- 
tional states of the products — are among such. These in- 
creasingly sophisticated experiments demand finer statis- 
tics on the QCT results for highly resolved, intrinsically 
quantum features are much more difficult to 'converge' 
using classical mechanics. The 'brute force' approach 
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is of course to increase the number of trajectories in 
the classical ensemble — and as their propagation is eas- 
ily task-parallelised, i.e. particularly suitable to solve in 
computer clusters — running millions of such is becoming 
common. 

In this work we report on an alternative QCT algo- 
rithm for investigating polyatomic processes which ef- 
fectively avoids any binning or weighting method while 
provides both vibrational and rotational resolution. Its 
main advantage comes from the fact that initial condi- 
tions are directly generated at the semiclassical rovibra- 
tional states of interest in terms of angle-action variables, 
which automatically reduces the total number of trajec- 
tories to integrate. These are however propagated in 
Cartesian coordinates, by means of a transformation very 
recently published In addition, a quasi-classical for- 
mula is derived to further incorporate possible rotational 
structures. To introduce and illustrate the capabilities 
of our method, we simulate highly-structured empirical 
correlated product distributions from the dissociation of 
ketene at 308 nm 

The paper is organized as follows: the problem, the- 
oretical aspects of the QCT implementation and results 
of its application are dealt with in Section |ll] and later 
discussed in Section [nil at last, Section Hv] summarizes. 



II. THEORY AND APPLICATION 



Ketene dissociation at the So electronic surface 



The photo-fragmentation of ketene. 
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present, it is firmly established that following optical ex- 
citation to the A^A" state, the molecule either undergoes 
intersystem crossing to the lowest-lying, triplet, or fast 
internal conversion to the singlet electronic state. From 
these, dissociation into methylene and carbon monoxide 
occurs. Despite the triplet threshold lies ~3150 cm~^ be- 
low the singlet, the fact that it presents a (small) barrier 
to dissociation — of a few cents of inverse centimetres — , 
but mainly, that the interaction which accounts for the 
change in electronic spin is rather weak, makes the sin- 
glet channel statistically dominant from excess energies 
as low as 100-200 cm^^. Such conditions make this pro- 
cess an effective prototype for a barrierless polyatomic 
unimolecular reaction on a .single potential energy surface 
(PES). Experimentally, it has the additional advantage 
of a well-defined total energy and time origin distinctive 
of photo-dissociations p!3 |. 



Earlier measurements of energy-dependent rates and 
overall product distributions were used as benchmarks 
to test and develop various statistical models including 
phase space theory (PST) and Rice-Ramsperger-Kassel- 
Marcus (RRKM) theory, cf. IH, [Ij, [23|. It was 
concluded that under moderate excitations the molecule 
dissociates statistically, exit-channel effects becoming in- 
creasingly important with augmenting energy. The dis- 
sociating picture therefore evolves from that of a loose to 
a tight transition state (TS) [3l|. During the last decade, 
however, there has been growing interest in the problem 
of correlated product distributions [l3, [13, l33 . In the 
most recent experiments, Komissarov et al. [1J| studied 
the photo-fragmentation of the molecule at 308 nm and 
measured product translational energy distributions for 
specific rotational states of CO, i.e. P{Et; jco)- These 
are the observables we will focus on in the present work. 



B. Singlet potential energy surface 

Extensive trajectory calculations require a reliable nu- 
merical representation of the singlet PES. This was gen- 
erated from the most recently reported high-level ab ini- 
tio data [1§|. In this work, Klippcnstein, East and Allen 
used a separation in transitional and conserved modes to 
determine the quadratic force fields along the reaction 
path — at the CCSD/[5s4p2d, 4s2p] level of theory — in 
terms of Cj symmetry internal coordinates. For excess 
energies of the order of that obtained with a 308 nm laser 
{Ecxc ~ 2350 cm^^), these force fields should provide all 
essential information on the couplings between the intra- 
fragment modes and the dissociation coordinate in the 
vicinity of the reaction path, as well as the interaction 
dependence on the mutual orientation of the fragments^. 
To further assess its quality, in the same work, the au- 
thors showed that variational RRKM energy-dependent 
rates predicted with this PES are in quantitative agree- 
ment with the experiment. 

Being S = {Si} the symmetry internal coordinates, 
Fij {i,j = 1, 9) the quadratic force constants and Vcc 
(Fig. 7 in Ref. [23|) the minimum potential energy along 
the reaction path as a function of the C-C interatomic 
distance, S4, the PES may be written as 



V{S)^Vcc{S4)- 
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where most Sok — SokiSi) can be straightforwardly cal- 
culated from the definitions of symmetry internal co- 
ordinates and the analytic parametrizations for the re- 
action path (Tables VI and III, ibid.). The slightly 



^ We discuss the implications of this quadratic approximation 
to the PES for trajectory calculations in some detail in Sec- 
tion [nni] 
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more elaborate are Sqq = \/2arccos [— cos (a/2) COS7] 
and 5*09 = arccos [— cos {a/2) sin7 sin~^ 0\ , where — 
if individual ketene atoms are numbered as they enter 
in the chemical formula — a — Z(H^C^H^), /? = /3i^2 = 
Z(H-^'^C^C^) and 7 (CH2 wagging angle) are all angular 
internal coordinates defined in Fig. 1 of the same work. 
The Fij functions were obtained from fittings to values 
in Table VI, ibid. 

Finally, we validated the PES by diagonalising the Hes- 
sian matrix and confirming that all local harmonic fre- 
quencies previously reported (Table VII in Ref. 29]) were 
successfully reproduced. 



C. Quasi-classical trajectories 

We will introduce in what follows a particular im- 
plementation of QCT ideas and show its predictive 
power in the simulation of highly-resolved jco-correlated 
translational energy distributions, from the unimolecu- 
lar dissociation of ketene. Theoretically, it can be as- 
sumed that under the experimental conditions reported 
in Ref. these product distributions correspond to a 
given total energy E = i^oxc + -E-nmodcs and angular mo- 
menta J, i.e. P{Et;jco) = P{Eu E, J, jco)- The zero- 
point energy in the molecule normal modes, -Enmodcs ~ 
4730 cm-i. 

In reading the following section, it is important to keep 
in mind that we use an association rather than a dissocia- 
tion picture for propagating trajectories. This is based on 
a method due to Hamilton and Brumer [35j which will be 
introduced in Section III C 21 where its many advantages 
are also discussed. 



1. Initial conditions: Angle-action variables and the 
transformation to Cartesian coordinates 

Within QCT, it is always preferred to generate the en- 
semble of classical initial conditions for the very precise 
quantum states explored in the experiment. Such a gen- 
eration is natural ly p erformed in terms of angle-action 
variables [H, 113, the actions being the classi- 

cal counterpart of quantum numbers^. These variables, 
however, should preferably not be used to run trajec- 
tories for they may lead to strong numerical instabili- 
ties. In fact, Cartesian coordinates are ideal to solve the 
classical equations of motion for molecular processes in- 
volving polyatomics. Codes built on these coordinates 
should result easily scalable and standard numerical ap- 
proaches, adequate. Therefore, the transformation from 
angle-action variables to Cartesian coordinates is an im- 
portant step in QCT. Very recently, we reported on such 
a transformation [T3| that fits to our purposes here. 



Discarding all center-of-mass coordinates the phase- 
space results 

r = {R,qi,a, (3,ai,ai,a2,ak,^i} 

(g){P,Xi,J,J^,l,ji,j2,k,Ki}; * = 1,4, (2) 

where R and P are the dissociation coordinate — the CH2 
and CO fragments centres-of-mass relative distance — and 
conjugate momentum. The remaining 22 angle-action 
variables are [13i | : qi, Xi, the vibrational phase and corre- 
sponding action of the ith normal mode of CH2 (i = 1, 3) 
and CO (i = 4); J, J^, the moduli of the total angular 
momentum and its z-axis projection in the laboratory 
frame; I, the modulus of the orbital angular momentum; 
ji, J2, the moduli of CH2 and CO rotational angular mo- 
menta; fc, the modulus of the total rotational angular 
momentum fc = Ji + J2! and ki, the algebraic value of 
the projection of j\ onto the C2 axis of CH2. Finally, a, 
(3, ai, ai, a2, ccfe and 71 are respectively all the angles 
conjugate to the angular momenta mentioned above. It 
is apparent that in terms of these variables we can effec- 
tively set up, individually, all the vibrational and rota- 
tional modes of the fragments. 

Given E, J (fixed to 1) and j2 (= Jco) — for each set of 
energetically permitted vibrational states Xi, i = 1,4 — 
the total available energy -Eavaii and maximum allowed 
values /max and ji^inax are completely determined'^. From 
all these, J^, k and ki are thus randomly selected from 
uniform distributions according to e [—J, J], |ji,max — 

j2\ < k < jl.max + j2 and Ki £ [- jl,inax, jl.max] • MicrO- 

canonical distributions are finally obtained by aleatory 
selection of all angles and phases from a uniform dis- 
tribution in [0,27r]. All vector constraints on the set of 
angular momenta must also be satisfied. 

At this point, the F vector is transformed to Cartesian 
coordinates, using the algorithm presented in Ref. p^ . 
Here, a slight readjustment is necessary due to the con- 
ditions of the experiment, as the total energy of the 
system — rather than the relative velocity between the 
fragments — must be fixed. Specifically, the three Carte- 
sian components of P (Section II C 3 in Ref. [3]) should 
be computed right before the final step (determination 
of nuclear momenta, second part of Section II C 11, 
ibid.) for only at this point it is possible to calculate 
the energy available to translation, -Eavaii.t, from which 
P — •\/2/^i?avaii,t (m being the system reduced mass). 
Computationally, an incubation period is needed before 
using an otherwise acceptable initial state to avoid pos- 
sible correlations in the numerical generation of pseudo- 
random numbers. 



^ Throughout this work all actions will be quoted in h units. 



^ The choice of ij will become apparent in Sections lII C 2l and lII C 3] 
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2. Propagation: Exit-channel corrected phase space theory 



As previously seen f Section III A|) . it is possible to as- 
sume that the photo-excited ketene forms a long-lived 
intermediate complex prior to fragmentation, which jus- 
tifies the use of a microcanonical distribution at the TS. 
In this work, extensive trajectory calculations were per- 
formed in conformity to the so-called exit-channel cor- 
rected phase space theory (ECCPST), proposed in the 
1980s by Hamilton and Brumer [35|. Within this method, 
a microcanonical set of trajectories are initiated in the 
products and directed towards the unimolecule region, 
the statistics being performed with those 'successful' tra- 
jectories reaching the TS. Clearly, the underlying as- 
sumption here is that the time-reversed dynamics of this 
subset of trajectories yields a microcanonical distribu- 
tion at the TS [IBl- The method bears certain similari- 
ties with the statistical quasi-classical trajectory (SQCT) 
model for atom-diatom insertion reactions of Aoiz and 
co-workers [ioj . 

Besides its simplicity, the ECCPST offers quite a few 
practical advantages, e.g. (1) initial conditions are gen- 
erated where the assumptions for partitioning into prod- 
uct normal modes are best fulfilled — i.e. at the inter- 
fragment distance for which the interaction between the 
CH2 and CO is negligible — so the usual dispersion in the 
ensemble's total energy will be minimum (of about 2- 
3 cm~^ in our calculations); (2) all threshold behaviour 
is expected to be reproduced at its best, since product 
states are given precise energies in direct correspondence 
to their quantum counterparts; (3) initial conditions are 
straightforwardly defined by fixing the centres-of-mass 
distance between the resulting fragments — which explic- 
itly appears in the transformation from angle-action vari- 
ables to Cartesian coordinates [l3| — while the TS can 
be re-defined at convenience during propagation (in this 
case, by the C-C interatomic distance which was proven 
best [23); (4) since product states statistics is performed 
on a subset of initial conditions, the extra (non-trivial) 
work of discretising classical magnitudes at the end of 
the trajectories is avoided; (5) any binning or weighting 
to cope with the conceivable non-adiabaticity of vibra- 
tions in the course of reaction — which was found to be 
of about 15-20% after some preliminary calculations — 
can be avoided; and (6) the PES does not need to be 
determined within the reactant region. A few other ad- 
vantages will be discussed in the following sections. 

The main drawback, however, is that (possibly) many 
trajectories are propagated which never actually reach 
the TS. In fact, we found that approximately only 1 out of 
10 trajectories enters the reactants region"'. The overall 



* The indetermination on the exact TS location makes necessary to 
discard some trajectories even when propagating outwards from 
the reactants region. In this case, these constitute nearly one 
half of the total ensemble. 
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FIG. 1: (Colour online) Vibrational branching ratio for singlet 
CH2 in coincidence with selected GO{x4, = 0, 72) states. Com- 
parison with experimental and previous PST calculations. 



propagation time may thus be significantly reduced by 
applying strategies to quickly identify those trajectories 
which will ultimately be discarded. 



3. Computational details: integration and convergence 

Except when stated contrary, batches of 5 x 10** tra- 
jectories were started at = 10 A and propagated 
towards the unimolecule region in accordance with EC- 
CPST f Section m C 2\ . Convergence of reported product 
translational energy distributions was tested by confirm- 
ing that those obtained with 2 x 10^ trajectories are sta- 
tistically identical. 

A fixed-step fourth-order Adams-Bashforth-Moulton 
predictor-corrector method initialised with a fourth-order 
Runge-Kutta integrator was used for propagation. Tra- 
jectories were followed until one of two conditions was 
fulfilled: (1) the molecule reached the transition state (at 
the C-C interatomic distance = Rcc, as relevant to 
the specific excitation energy from Table X in Ref. ^29]); 
and (2) t = T^ax = 20 ps. The step size At = 2x 10^^ ps 
was adjusted so that the maximum error in total energy 
did not exceed AE = 20 cm"'^, i.e. less than 0.3% of E, 
with which the mean error was (AE) = 5 cm^^. Under 
these conditions, the fraction of non-fragmented to frag- 
mented trajectories was consistently smaller than 0.5%. 
Typically, a single trajectory required between 1 and 3 
seconds of CPU time — according to the specific values of 
X2 and j2 — on an Altix XE210 computer with two Dual- 
Core Intel® Xeon® 5030 ('Dempsey') processors. 



4. Correlated translational energy distributions 

After the subset of trajectories reaching the TS is iden- 
tified, the P'~' {Et; E, J,X2,jco) distribution correspond- 
ing to each X2 scissor mode of CII2 is calculated by stan- 
dard statistics over it. From these, the overall corre- 
lated translational energy distribution P'~'{Et;E, J, jco) 
is obtained by merging the different X2 contributions al- 
together, with the corresponding weights according to 
the vibrational branching ratios shown in Fig. [TJ 
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FIG. 2: (Colour online) Correlated product translational energy distributions after excitation with a 308 nm laser: QCT results 
with vibrational resolution. Comparison with the experiment. 



In principle, to predict these weights at a given j2 state 
of interest two contributions need to be taken into ac- 
count. First, Monte Carlo methods should be applied to 
find the relative probability that for each scissor state X2 
a given set of angle-action variables is finally accepted 
as a valid initial condition. Here the ECCPST proves 
again particularly convenient, since this is readily done 
if track is kept on the percent of successful trials and 
the hypervolume in phase space from which each X2 set 
of initial conditions was randomly selected. The second 
step is to correct this relative probability for dynami- 
cal effects, which in practice would mean to determine 
the ratio of capture probabilities from each scissor mode. 
Despite its apparent simplicity, obtaining a reliable ap- 
proximation for the latter could be difficult for it not only 
depends on the quality of the PES, but also on the rel- 
ative accuracy of TS locations and the set of integration 
parameters. However, it was previously shown that, by 
cancellation of errors, least-biased values (also referred 
as PST, depending on the convention em ploy ed) are in 
excellent agreement with the experiment [ij]. The val- 
ues reported in Fig. [T] are directly calculated from our 
initial conditions by the algorithm explained at the top 
of this paragraph. Their agreement with experimental 



and previous theoretical values confirms the correctness 
of our approach. At last, the theoretical distributions are 
Gaussian-blurred to match the experimental resolution. 

Calculated curves for several representative j2 states 
are compared with measurements in Fig. [2] — and have 
been quoted as 'A = 0.0' (blue), the reason for which will 
be given in the following section. The reported contri- 
butions from adjacent rotational states |14l | was appro- 
priately taken into account to produce the distributions 
corresponding to j2 = 4 and 9. Additionally, the zero- 
point-energy (ZPE) corrected distribution [41] has been 
included in the j2 = 4 panel, which illustrates the gen- 
eral behaviour for all other rotational states. To produce 
the latter, an ensemble of 2 x 10^ initial conditions with 
continuous actions was started in the products, being 
the vibrational modes constraint between zero and their 
maximum energetically allowed values. 



5. A modification to the PES 

Except of course for the j2 — 24 rotational state, all 
'A = 0.0' theoretical curves in Fig. [5] have the correct 
bimodal shape due to the excitation of the first and 
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the ground states of methylene's scissor mode (with fre- 
quency W2 ~ 1400 cm"'^). The experimental positions of 
the peaks are also fairly well reproduced. It is important 
to stress here that neither of these facts is trivial. Note, 
for example, the strong discrepancy of the ZPE-corrected 
curve which fails to reproduce any of them. Moreover, to 
achieve a similar agreement with the experiment in the 
case of the triatomic dissociation NCO — > N + CO, the 
GM method had to be applied [9[. However, the diver- 
gence between the experimental and ZPE distributions 
has an additional important implication. As the popula- 
tion corresponding to the ground state of the CH2 scissor 
mode 0:2 = is practically negligible, even by using the 
GW one could never obtain the proper translational en- 
ergy distributions from ZPE-corrected calculations only. 
It is thus crucial that, at this point, the 'vibrational' res- 
olution has been naturally obtained within QCT and no 
binning or weighting procedure is needed. 

However, it is clear that within each vibrational lobe, 
these theoretical distributions tend to predict larger than 
observed populations of low CH2 rotational levels. A 
plausible explanation may be derived from calculations 
on triatomic complexes. From these, it is known that 
quadratic representations of PES systematically tend to 
overestimate the energy transfers from the diatomic rota- 
tion to the inter-fragment translational degree of freedom 
[42] . Physically, this is due to the bending mode being 
unrealistically hindered, i.e. insufficiently damped, in the 
course of fragmentation. 

The previous reaction rate calculations on this PES 
[2^ , being an observable of 'global' character, would very 
likely not depend as strongly as the distributions we in- 
tent to reproduce. Nevertheless, only the construction 
of a new global PES or ab initio dynamics calculations 
could give a reliable answer to this issue. Both these are 
beyond the scope of the present work. 

On the other hand, several modifications can be at- 
tempted in order to bring the rotation-to-translation 
transfers closer to the empirical observations, under vari- 
ous physical grounds. One possible direction is to slightly 
change the force constants related with the molecular 
modes directly responsible for such transfers. The alter- 
native followed here was to rewrite the PES as 

Vx{S)=xV{S) + {l-x)V°^{S), (3) 

where x = exp [-A(5'4 - R^)] and V°° = Vcn^ + ^co is 
simply the asymptotic potential interaction. Eq. ([3]) con- 
serves both the PES at the TS and the correct asymp- 
totic behaviour. In particular, A = corresponds to the 
original PES. Physically sensible values for the damping 
parameter A, from calculations on triatomic systems, are 
usually between 2 and 3. 

Trajectories were propagated on the modified PES 
of Eq. ^ and the theoretical curves corresponding to 
A — 0.75 and 1.5 are included in Fig. [2l From these, it 
is clear that reducing the energy transfer between CH2 
rotation and inter-fragment recoil translation invariably 
leads to the overall experimental trend. The change is 



more strongly pronounced at higher CO rotational exci- 
tations and all A = 1.5 curves are in very good agreement 
with the experiment. 



D. 'Inserting' the CH2 rotational structure 

From Fig. [2] is readily seen that the main qualitative 
difference between ours {e.g., A = 1.5) and the exper- 
imental results is the rotational structure observed in 
the latter, which relates to that of CH2(a^^i). In this 
section we show how these features can also be theoret- 
ically reproduced by the adequate modification of the 
P^{Et;E,J,jco) distributions. 

1. A quasi-classical formula to get the rotational resolution 

Let US first consider the simpler though instructive case 
of a triatomic unimolecular dissociation ABC — > A -|- 
BC(u, j), at fixed total energy E and angular momentum 
J. Here, w,j represent the vibrational and rotational 
states of the diatomic fragment. In particular, we are 
interested in the vibrationally-resolved product transla- 
tional energy distribution (at a sufficiently large value of 
the inter- fragment distance 

Classically, after removing all center-of-mass coordi- 
nates, this is exactly 

5 

p^{Et;E,j,v)^ / nA,x(r)dr, (4) 

being 

r = {R,q,a,P,ai,a-i,P,x,J',J^,l,j} (5) 

the physically relevant set of angle-action variables (com- 
plemented with the dissociation coordinate R and its con- 
jugate momentum P, e.g. Eq. ([2])), and 

Ai S{nt-Et), (6) 

A2 = d{n-E), (7) 

A3 = S{J'-J), (8) 
A4 = 6{x-v), (9) 

A5 = -Q{P)6{R-R^), (10) 

where Ai to A4 explicitly account for all fixed magni- 
tudes while A5 limits integration to the outgoing flux 
through the R = R°° surface. As usual, S and O are 
respectively Dirac's 5- and Heaviside functions. In addi- 
tion, the total Hamilton function has been partitioned as 
Ti. — Ti.t+Ti.„, i.e. into its translational and rovibrational 
components, with 

H. = i^(p^ + ^)+y.,BC. (11) 
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By definition, the interaction potential of A with respect 
to BC, Va,bC: tends to zero as i? ^ R°° . 

At last, x(r)dr is the probability for the system me- 
chanical state to be in the volume T of phase space. The 
x(r) function can be rewritten as 



x(r) = Cjjj.xir), 



(12) 



where Cjjjj — A{J,l,j)Q{J ~ \ Jz\) holds all vector 
constraints, A representing a triangular inequality, x 
is effectively independent of (provided this is large 
enough) and carries all the information on the full disso- 
ciation process: that coming from both the initial state 
distribution and possible dynamical effects, i. e. the state- 
to-state transition probabilities. 

On recalling some basic (5-function properties, partial 
integration in Eq. ^ can be shown to yield (we sim- 
plify notation by discarding the explicit parametric de- 
pendence with respect to E and J in what follows) 

P'^iEuv) cx J Cjjj,6 {E-E,- E„) 

xx{Et,J,,l,j;v)dJAldj, (13) 

where Ej-^ is the rovibrational energy of BC and x is the 
angle/phase average of x- The centrifugal orbital term 
has been neglected for it depends on which becomes 
negligible as i? — > R°°. Using the mean value theorem, 
this can be expressed as 



P'^{Et;v) = E'^iEt;v)Pi:B{Et;v), 



(14) 



being 

PEsiEuv) = J CjjjjSiE -Et~ E,,)d.hdldj (15) 

the (classical) least-biased vibrationally-resolved product 
translational energy distribution. 

Eq. (|14|) may even look deceptively simple, though its 
derivation here proves useful to disentangle the physical 
origin of its factors. In a sense, the underlying ideas can 
be related to the theory of surprisal analysis of Levine 
and Bernstein [l^, if the least-biased plays the role of 
the prior distribution. Eq. p4| basically states that P'~^ 
differs by a scaling factor (which in turn depends on Et 
and v) with respect to the least-biased value P^g- The 
factor itself, S*^, 'inherits' from x the dynamical in- 
formation on the particular process, but it is additionally 
'weighted' by the density of constraint rotational states. 

In obtaining a quasi-classical analogue to Eq. two 
approximations are made. First, we replace the integrals 
in Eq. (jlSp by the corresponding quantum sums and sec- 
ond, we postulate that if QCT founding ideas are valid for 
the process under study, using 'EF in the place of its quan- 
tum analogue 2'3 must be a fairly good approximation. 
We should note that the latter is not only a dynamical 
statement, but also a subtle assumption on the equiva- 
lence between the classical and convoluted-quantum den- 
sity of rotational states, which is expected to hold under 



quite general conditions^. By means of these approxima- 
tions one can write 

P«^(£;t;«) = E''iEuv)P^siEuv), (16) 

which explicitly reads 

(17) 

after summing over all possible final vibrational states. 

Although our derivation has been made for a triatomic 
system, it is easy to write the equivalent of Eq. pT|) 
for the correlated translational energy distribution in the 
fragmentation of ketene 



P'^''(Et;j2) 



E 



pC{Et;v,,j2) 



a, 



Vvi,J^,Lk,ji,T ' LB 

xS{E~Et~E„), 



(18) 



where 



Cjjjknj2r = A{J,l,k)A{k,ji,j2)e{J-\J,\) 

xe(ji-|T|). (19) 

that finally allows to insert the CH2 rotational structure 
to our vibrationally-resolved distributions in Fig. [2l Here 
ji and T constitute mere labels for the CII2 rotational 
energy levels, as will become apparent in Section Hi D 21 
The ECCPST proves again particularly convenient as 
the P^g distribution in Eq. can be directly obtained 
from the ensemble of microcanonical initial conditions. 



2. CH2(a} Ai) asymmetric rotor energy levels 

Methylene has different its three principal moments of 
inertia {i.e. laa < hb < Ice), being thus an asymmetric 
rotor (AR). There is no general quantum formula which 
determines the rotational energy levels for such molecules 
in terms of the set of (quantum numbers) action vari- 
ables, ji, Kj, used for computing the corresponding initial 
conditions [13 |. 

However, it can be proved that the rigid-rotor Hamil- 
tonian Tirot becomes block-diagonal in the matrix repre- 
sentation of symmetric rotor (SR) eigenfunctions [4^14^. 
It is also possible to show that for a given angular mo- 
mentum J, every AR state correlates with one prolate 
{J,Ka) and one oblate {J,Kc) SR level only. It is thus 
a very common practice to label the AR levels as JxaK^, 
the indexes taking the values Ka G {0, 1, 1, 2, 2, ... , J, J} 
and Kc € {J, J, ...,2,2,1,1,0} with increasing energy of 
the terms. As it is easily seen, one could alternatively 



This is not generally true for vibrational states, hence we deal 
with vibrationally-resolved distributions. 
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FIG. 3: (Colour online) Correlated product translational en- 
ergy distributions after excitation with a 308 nm laser: QCT 
results with rovibrational resolution. Comparison with the 
experiment. 



use a single label defined hy t — Ka — Kc, which in turn 
spans the 2 J + 1 values from — J to J. r rather than ki 
has been used in Eq. to stress the fact that this label 
has no direct correspondence to the ni action. 

In our calculations, spectroscopic CH2(a^Ai) energy 
levels were used when available and complemented with 
theoretical AR calculations at higher energies [i^, [47l |. 
The error in theoretical values could be of several inverse 
centimetres, especially at large values of ji and r. 



III. DISCUSSION 

Rotational structures are barely visible in the measure- 
ments for j2 = 4 and 9, being strongly blurred due to the 
contribution of adjacent rotational states. Fig. [3] then 
shows the application of Eq. (fT8|) to A = 1.5 theoretical 
distributions, corresponding to j2 = 15 and 24, as a test 
on the capabilities of the method. Within each panel, 
two theoretical curves correspond to convolutions with 
different strengths. One of them (darker, red) has been 
chosen to closely match the experimental resolution [l3| 
while the other (blue) shows the typical fine structures 



which are attainable with our method. 

The empirical rotational structure is fairly well repro- 
duced by the QCT calculations, being all theoretical and 
experimental peaks correcly aligned. The larger disagree- 
ment is observed within the upper panel, in the region 
shown by an arrow. This area corresponds to the match- 
ing between the low energy tail of the X2 — state and 
the high energy component of X2 = 1, being the former 
the main contribution to it. The discrepancy between our 
curve and the experiment is thus a secondary effect due 
to the overestimation of rotation-to-translation transfers 
already discussed in Section fll C 51 It cannot be seen in 
Fig. [5] as the convolution used there artificially 'popu- 
lates' the missing rotational states. 

The least-biased distribution corresponding to j2 — 24 
is also depicted in Fig. O obtained from the whole set of 
initial conditions. This distribution is directly compara- 
ble with the PST calculations shown in Fig. 8 of Ref. p^ . 
In fact, application of Eq. (|18p in order to 'add' a ro- 
tational structure to this curve would trivially reduce it 
to the quantum version of the least-biased/PST distribu- 
tion. At last, by comparing this with the final theoretical 
curves, it is readily seen how the suppression observed in 
the low-i?t part are correctly described by QCT. 

It is important to address one last point referring to 
the insertion of the rotational structures in Section IIIDI 
In principle, one might think that a different, simpler ap- 
proach could be used: To additionally fix all rotational 
quantum numbers and propagate separate ensembles for 
each semiclassical state, very much as we have done here 
for each scissor excitation or can be found elsewhere, 
e.g. po| . However, with respect to the initial conditions 
there is a basic difference in the way the vibrational and 
rotational modes are determined [1J| . Using the normal- 
mode approximation, the former are assigned precise en- 
ergies corresponding to the quantum values. This means 
that normal mode displacements and conjugate momenta 
are calculated to exactly match the desired energy in the 
particular mode. On the other hand, classically, after fix- 
ing the angular momentum vector for a particular frag- 
ment, there remains an indetermination in the rotational 
energy due to the vibration-rotation coupling. In other 
words, there are different inertia tensors due to dissim- 
ilar normal-mode displacements which will result in dif- 
ferent energies even for the same angular momentum. In 
the case of ketene, for example, calculations after addi- 
tionally fixing the values of ji and ki yield broad struc- 
tures in the correlated translational energy distributions 
which are several hundred inverse centimetres wide. This 
would completely blur the rotational structures observed 
in Fig. [3) The indetermination should be smaller for sym- 
metric tops — where a formula exists that determines the 
energy in terms of the angular momentum and one of its 
components — but still would, in general, strongly limit 
the QCT rotational resolution. 
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IV. SUMMARY AND CONCLUSIONS 

A QCT algorithm has been presented that aUows to 
simulate indirect polyatomic processes involving exit- 
channel effects up to rovibrational resolution without the 
use of any binning or weighting procedure. It consists of 
three stages: (1) initial conditions are started in terms 
of angle-action variables to individually select the ex- 
perimental states to be studied, which are later trans- 
formed into Cartesian coordinates for trajectory integra- 
tion; (2) trajectories are propagated using an inverse dy- 
namics simulation in the spirit of the exit-channel cor- 
rected phase space theory of Hamilton and Brumer; and 
(3) finally, an approximate quasi-classical formula is used 
to accurately incorporate possible rotational structures 
into the quasi-classical product distributions. Of these, 
only the second step is linked to the statistical assump- 
tions inherent to the study of indirect processes, being 
the other two of general applicability. 

The method was introduced and its capabilities il- 



lustrated with the simulation of highly-resolved corre- 
lated translational energy distributions from recent ex- 
periments in the photo-dissociation of ketene at 308 nm 
[14| . After slightly modifying the PES to deal with an 
overestimation of the energy transfer among rotational 
and translational degrees of freedom, calculated distri- 
butions with both vibrational and rotational resolution 
are in good agreement with the experiment. 
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